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The hydrophobic effect (HE) is commonly associated with the demixing of oil 
and water at ambient conditions and plays the leading role in determining the 
structure and stability of biomolecular assembly in aqueous solutions. On the 
molecular scale HE has an entropic origin. It is believed that hydrophobic 
particles induce order in the surrounding water by reducing the volume of con- 
figuration space available for hydrogen bonding. 

Here we show with computer simulation results that this traditional picture 
is not correct. Analyzing collective fiuctuations in water clusters we are able to 
provide a fundamentally new picture of HE based on pronounced many-body 
correlations aff"ecting the switching of hydrogen bonds between molecules. These 
correlations emerge as a non-local compensation of reduced fiuctuations of local 
electrostatic fields in the presence of an apolar solute. 

The HE has a multifaceted nature, i. e. its physical manifestation depends on the length- 
scale^. On the mesoscale, i.e. hydration of an assembly of hydrophobic units or an ex- 
tended hydrophobic surface, HE is driven by energy/enthalpy and occurs as a 'dewetting' 
transition-^^ which has far-reaching consequences for processes such as protein folding-^iS and 
nanoparticle self-assembly-. Meanwhile, HE on the molecular scale has an entropic origin^i^, 
particularly near room temperature and lower, while it is believed to eventually become 
energy /enthalpy driven at higher temperatures^»i^. Furthermore, the molecular scale hy- 
dration thermodynamics (hydrophobe solubilities, partitioning of hydration free energy into 
energy/enthalpy and entropy contributions, etc.) appear to be well established and can be 
worked out, for example, using scaled-particle theory^: or information theory^^. While those 
theories are successful in predicting solubilities of hydrophobic solutes and several related 
thermodynamic features, they do not provide deeper insight into the physical mechanism 
underlying HE. Notwithstanding all efforts and advances in the fielcii^-i^iiiii^^— the physical 
picture of HE is still far from being complete and even fundamental issues, such as the 
mechanism underlying hydrophobicity on different length scales, still have to be clarified. 
From the physical point of view the most puzzling feature of HE remains the microscopic 
picture of entropy loss upon hydrophobic hydration. 

Intuitively, the entropy loss is usually attributed to the reduction of configuration space 
available for hydrogen-bondingi*^, which is due to the fact that water molecules need to re- 
organize around a hydrophobic solute to avoid sacrificing hydrogen bonds. This is supposed 
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to lead to remnants of clathrate structures^, which are, however, not rigid and their quan- 
titative importance for understanding hydrophobicity remains questionable^^. Furthermore, 
even the actual role of hydrogen bonds for the HE is apparently not entirely clear— li^. 

Thus, the physically most intriguing question to be answered still remains: If HE is 
entropy driven, what specifically causes the loss of entropy? Since entropy loss is directly 
related to a reduction of available volume in configuration space, how does it affect degrees 
of freedom of water molecules? 

Here we present conclusive simulation results which unravel a fundamentally new pic- 
ture of the mechanism of HE based on pronounced many-body correlations affecting the 
intermolecular exchanging of hydrogen bonds. We carry out constant pressure Monte-Carlo 
simulations of TIP5P water^^ and model hard sphere solutes in an orthogonal simulation 
box with periodic boundary conditions (see Methods section for details). 

In order to present the conceptual change in our understanding of the HE we first 
address the inconsistencies of the traditional picture. The radial correlation function, 
= g(r) = p^{Yl!i=i ^l'" " kil))' is used to quantify the degree of translational ordering 
of water molecules around hydrophobic solutes. We find a non-monotonic dependence of 
the contact density on particle size, which is due to commensurability of the solute surface 
and water packing. The relative density fluctuations monotonically decrease with solute 
size. Two water molecules are defined to be in close contact if their intermolecular distance 
is less than 3 A , and they are said to be hydrogen-bonded if they are in close contact 
and if the angle O — H ■ ■ - O is larger than 150°. The cutoff distance for neighbors in close 
contact is set at 3 A and is more appropriate with respect to the conventional definition of 
3.5 A^^, as there is no preferential mutual orientation beyond the distance of 3 A (see Fig. 

The distribution of the number of water molecules in close contact and the number of 
hydrogen bonded contacts per water molecule located in the first and second hydration shell 
and in bulk water is shown in Fig. [T]d and c. Except for the smallest solute (with radius 
1.4 A) there is no appreciable difference (say of the order of > 0.5) in the number of total 
and hydrogen-bonded contacts with respect to bulk water, neither in the first nor in the 
second hydration shell. If the main effect of a hydrophobic solute would be the reduction 
of the configuration space for hydrogen bonding, then one would naturally expect to find 
less neighbors in close contact. Clearly, this is not the case. Moreover, the distribution is 
much narrower in bulk water, which already suggests that the small-scale fluctuations in the 
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vicinity of hydrophobes (i.e. hbrations, hydrogen-bond exchange, etc.) are enhanced with 
respect to the bulk. According to our criterion for the nearest neighbor we find that the 
probabihty of a water molecule having 3 hydrogen-bonded nearest neighbors is negligible, 
irrespective of its position. 

Aside from an altered number of nearest neighbors the traditional picture also suggest 
a more ordered local structure. In order to cause entropy loss the structural fluctuations 
should tend to diminish. To inspect in detail the structural ordering of water molecules 
in close contact, around hydrophobic particles and in bulk water, we employ the recently 
introduced dipolar order parameter—. 



where i ^ j and the sum is taken over the neighbors in the first coordination shell of the 
i-th water molecule and a^"^ and correspond to maximal and minimal dipole-dipole 
potential at a given intermolecular unit vector and unit dipole vector of the tagged molecule. 
This way, D{i) takes values between (maximal repulsion) and 1 (maximal attraction). The 
distributions of dipolar ordering shown in Fig. |2] confirm the idea that water is on average 
slightly more orientationally ordered in the first hydration shell of hydrophobic particles 
(note the slight shift of the distributions towards higher D values). More importantly, the 
width of the distributions p{D), which reflects the constraining of orientational degrees of 
freedom, remains rather unaffected. This width clearly increases upon increasing polarity 
of the solute indicating an orientational relaxation in the hydrogen-bond network—. The 
situation in the second shell is similar, albeit less pronounced. Most likely, due to a commen- 
surability effect, the second shell molecules exhibit a reverse trend in orientational ordering 
in the case of s2. 

Apparently the pair interactions and consequently also the hydrogen bonds are strength- 
ened but there is no reduction of the orientational configuration space explored by individual 
water molecules in the hydration shells as the distributions are merely shifted while their 
form remains unchanged. Thus, we find that i) there are no significant differences in the 
number of nearest HB and non-HB neighbors and ii) there is no orientational constraining 
although the local structure is more ordered. This is clearly in contradiction with the idea of 
reduced configuration space available for hydrogen bondingi*^. So how can this contradiction 
be reconciled? 
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There is an important difference between the total volume of configuration space and the 
volume that is actually visited at a given temperature. As the number of close contacts 
is not altered significantly by the mere presence of the solute the latter is not expected to 
significantly affect the volume of the configuration space available to a water molecule for 
hydrogen bonding. On the other hand, in the liquid state water molecules form a dynamic, 
labile hydrogen-bond network. There is clear evidence that the molecular reorientation 
underlying the exchange of hydrogen bonds is not diffusive {i.e. does not occur as angular 
Brownian motion) but rather proceeds in terms of large-amplitude sudden jumpsi^i^"— . As 
it is impossible for water molecules to each form four HBs at a time in the liquid state, 
this introduces a strong frustration into the system emerging from the competition between 
the entropic and energetic driving forces. This means, that while fluctuating about the 
minimum energy conflguration exhibiting e.g. librational movements, a water molecule is 
most likely to form two hydrogen-bonds to its closest neighbors. As a water molecule in the 
bulk liquid is intrinsically inclined to form more than one hydrogen bond, this necessarily 
leads to conflnement of its orienational conflguration space. This fact suggests to rationalize 
the entropy loss in HE in an alternative way. 

The current understanding of hydrogen-bond network dynamics in neat water already 
anticipates the collective concerted motions of several water molecules^. This suggests 
looking for the potential source of entropy loss in the altered 'communication' between 
emerging transient water clusters involved in the exchange of hydrogen bonds, i. e. the re- 
duction of conflguration space due to constraining of collective degrees of freedom. For large 
enough fluctuations the transient water cluster can rearrange into another conflguration just 
by the exchange of hydrogen bonds. Hydrogen-bond fluctuations (for example, librations) 
try to compensate for the entropic frustrations caused by a transient mutual alignment of 
water molecules in between sudden jumps. Since structural fluctuations are caused by lo- 
cal fleld fluctuations they are expected to be dramatically suppressed in direct vicinity of 
a hydrophobic solute, because the latter exerts only an extremely weak electrostatic fleld. 
Without an additional compensation mechanism the hydrogen bonds would be signiflcantly 
strengthened and the entropic frustrations would be expected to grow further. Therefore 
there must exist a tendency of nearby water molecules (nearest and next-nearest neighbors) 
to compensate for the suppressed fluctuations of local electrostatic flelds. This in turn can 
not happen unless the fluctuations of nearby hydrogen-bonded clusters (intra- and interclus- 
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ter fluctuations) become correlated as to maximize local field fluctuations while maintaining 
as many as possible mutual water arrangements close to the optimal HB geometry. Thus, in 
order to satisfy the local energetic demand to form hydrogen bonds, the resulting entropic 
frustration relaxes non-locally in the vicinity of a hydrophobic solute. 

Testing this hypothesis demands the evaluation of various many-body correlations. Al- 
though addressing the problem of coupled translational-orientational multi-body correlations 
is nontrivial it can be signiflcantly simplifled in the following manner. While single water 
molecules are, by nature, indistinguishable they can be transiently classifled as being hy- 
drogen bonded or non-bonded to its nearest neighbor (using the same criteria as above). 
Any HB exchange event merely permutes the indices between water molecules. With these 
criteria we can, at any instant, classify the HB and non-HB neighbors of any given molecule. 
Involving various types of fluctuations on different time-scales, such as librations, dimer 
tumbling and HB jumping (see for example^^ and references therein), following both the 
relative positions and orientations of water molecules individually is simply too complex 
and in fact unnecessary. To determine how molecular degrees of freedom are constrained in 
the liquid state we shall be primarily interested in structural fluctuations in water clusters. 
Such fluctuations can be easily described in terms of fluctuations of interaction energies of 
hydrogen bonded and non-hydrogen bonded molecules. 

Using the random variable transformation theorem^I we can map the joint conflgurational 
probability density onto its functionally dependent joint probability density for pair potential 
energies. Thereby a given matrix of molecular positions and orientations, f2, is mapped onto 
a vector of interaction energies, U, which can be formally written for given values uj and u 
as follows: 



The functional relation 6{u — f{u:^)) contains a class indicator (close contact, hydrogen- 
bonded, etc.) as well as the appropriate averaging operation accounting for the indistin- 
guishability of pairs within a certain class (for details see section 1 in the Supplementary 
information). This way we can construct, using appropriate functional relations, probability 
densities for observing, for example, a (non)hydrogen bonded pair with certain energy or 
a joint probability density of observing two pairs of different classes having given energies. 
Comparing the values per water molecule/cluster inside the flrst and second hydration shell 
with the corresponding bulk values we are able to select the most important contributions 
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to the entropy loss. The Gibbs-Shannon entropy, S[p] = —ks f p{ui) In p{ui)dui, is used to 
evaluate the total uncertainty of a quantity Ui, and thus quantifies the fiuctuations. Al- 
ternatively, Kullback-Liebler entropy or Correlation entropy (the latter is taken after--) is 
used to quantify the total correlation between two random variables X and Y (which can be 
components of U, for example various combinations of hydrogen-bonded and non-HB pairs 
{Uhb and Unon-HB)- The correlation entropy can be expressed in terms of individual and 
joint entropies, Sc{X, Y) = S{X) + S(Y) — S{X, Y), or exphcitly in terms of corresponding 
probability densities: 



We limit the present discussion to 3 and 4-particle correlations to asses how fiuctuations of 
various HB and non-HB pairs are correlated. 

Being predominantly interested in generic features and less on the specific effect of solute 
sizes we focus first on the properties of individual pairs of molecules and find that the dipolar 
entropy difference between the hydration shells and the bulk, which directly measures the 
degree of orientational constraining, does not show a monotonic behavior with respect to 
the solute size and is less than 1% of the corresponding bulk value (for details see section 2 
in the Supplementary information). This confirms the assumptions that a reduced orienta- 
tional configuration space of single molecules is not the reason for HE. This fact is further 
substantiated with the values of the entropy difference, AStot-c, of the fluctuations of the 
total interaction energy of tagged molecules with its nearest neighbors, S[p{J2iUio)], (Fig 
[HH). It may also be taken as a contribution to the fluctuations of the local electrostatic 
fleld (solely) due to the thermal motion of nearest neighbors. In fact, AStot~c is mostly posi- 
tive, except for the smallest solute which indicates that local fleld fluctuations due to nearest 
neighbors are enhanced in the vicinity of the hydrophobe. Meanwhile, the entropy of mutual 
fluctuations of the potential energy of both hydrogen bonded neighbor pairs (Fig. |3^ blue) 
and non-hydrogen bonded pairs (Fig. [3^ green) decrease with respect to bulk water. The 
former can be understood as a measure of librational-type of fluctuations and its lowering is 
indicative of a slight HB strengthening^^ (though the relative difference is smaller than 2% 
of the bulk value). Together with decreasing fluctuations of the non-HB interaction energy 
and the simultaneous increase of fluctuations of the total interaction energy with nearest 
neighbors this immediately hints at higher order (beyond pair) correlations, as suggested 
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by our hypothesis. The local strengthening of HB is also suggested by the shift of p{D) 
towards larger D. The fact that fluctuations of the total interaction energy are so small and 
their entropy difference mostly even positive along with the strengthening of HB are fully 
consistent with the large- amplitude sudden jump picture of HB-network dynamics. A closer 
look at entropies of various 3- and 4-body correlations in Fig (Hb-d reveals striking, upto 22 
fold increases with respect to bulk water. The largest increase is observed in the case of 
correlations of nearby non-HB pairs bridged via a HB (12-22 fold), a HB pair and a nearby 
non-HB pair (7-9 fold 9-15 in the case of geminal and 9-15 fold in the case of vicinal pairs) 
and two HB pairs that do not share a common HB bond (are not bridged via a HB; 9-11 
fold). A significant, albeit smaller, increase is also observed in the case of geminal HB pairs 
(app. 100% increase) and HB pairs which are bridged via a HB (70-90% increase). Let us 
briefly mention that all 3- and 4-body correlations are, to linear order, anti-correlated. The 
pronounced many-body correlations propagate into the second hydration shell (for details 
see Table 1 in section 3 of the Supplementary information) but fall off rapidly further away. 
Thus, the picture of two perturbed hydration layers is retained. 

The results suggest that the traditional explanation of the HE on the molecular scale 
needs a substantial revision. It is not a reduced configuration space for hydrogen bonding 
that is responsible for the observed lowering of entropy, but a striking increase of many-body 
correlations (Fig. H^), which is essentially due to hydrogen bonding being a strong and ori- 
entationally dependent interaction. The increase of many-body correlations is necessary to 
compensate for the reduction of fluctuations in the local electrostatic field and the result- 
ing local HB-strengthening when one or more "polar" water molecules are replaced by an 
"apolar" hydrophobic particle. This effect should scale as the difference between local field 
distributions experienced by the water molecule in bulk (Fig. Hb black curve) and by a bulk 
water molecule for which one neighboring water molecule is omitted from the calculation of 
the local field (Fig. |Dd black curve). If left uncompensated, such as in the latter case, the 
distribution is shifted significantly to lower values. It also turns out that the distributions 
of field strengths in the hydration shells are almost identical to the one in bulk liquid (Fig. 
Hb symbols). Furthermore, by looking at the distribution of orientations of the local field 
over the unit sphere (for details see Fig. Ht and section 2 of the Supplementary information) 
shown in Figure HJi we find that the distribution is isotropic in the bulk as well as in the 
second hydration shell (not shown). Meanwhile, in the first hydration shell the local field 
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fluctuations are constrained to the plane containing the normal to the solute surface. Such 
a distribution ensures that all water molecules in the first hydration shell experience the 
maximum span of local fields and thereby local torques which counteract the entropic frus- 
trations caused by a transient mutual alignment of water molecules in between HB exchange 
events. Such a compensation mechanism comes at cost of many-body correlations. The local 
increased entropic frustration is thus compensated non-locally. 

The proposed picture of hydrophobic hydration does not imply any water immobilization 
in the spirit of the iceberg hypothesis. In explaining the moderate retardation of reorienta- 
tion dynamics2^i^ the current picture, based on the excluded volume effect, proposes that 
the slowdown is a result of fewer accessible configurations of the transition-state (TS) in 
a HB exchange event due to the presence of the solute. We demonstrate that the lower 
number of accessible TS configurations rather results from many-body correlations {due to 
correlations the TS is simply statistically less probable) rather than actual sterical hindrance 
of reorientation. Thus the expression for the slowdown factor remains unaffected. Our view 
also provides an explanation of the slightly better solubility of nonpolar solutes in heavy 
relative to light water'^^. While being more inert heavy water responds less strongly to lo- 
cal field fluctuations, which leads to HB strengthening and thus a higher cohesive energy 
density— but also results in a more structured liquid with respect to light water^^. Being 
intrinsically more structured the relative extent of local field fluctuations to be compensated 
upon the introduction of a hydrophobic molecule is therefore also smaller. 

Finally, the results allow us to speculate that the range of hydrophobic interactions 
beyond that expected from the minimally-exposed surface area reasoning is a result of the 
propagation of many-body correlations beyond the first hydration shell. 

I. METHODS 

We performed extensive constant pressure Monte-Carlo (MC) simulations with TIP5P 
water and freely moving hard sphere solutes at room temperature (298 K) and 1 atm assum- 
ing periodic boundary conditions. After an extensive equilibration period we performed as 
many trial moves as to assure that on average each molecule was successfully moved 4 x 10^ 
times at the maintained acceptance rate of 30%. To exclude (auto) correlations the successive 
configurations used for the analysis were taken to be as far apart as to assure that between 
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each taken configuration each molecule was successfully moved at least 5 times. The hard 
sphere radius of the water molecule when in close contact with the hydrophobic spheres was 
taken to be 1.4 A. The following radii of solutes were considered: r(sl)=1.4 A, r(s2)=2.1 
(1.5 X r^^^o) A, r(s3)=2.52 (1.8 x r^^^^) Aand r(s4)=2.8 (2 x r^^^^) A. According to the solute 
size the number of water molecules was (in ascending size) 2668, 2670, 3513 and 3560. The 
total number of MC steps was 35.57 x 10^ 35.60 x 10^ 46.8 x 10^ and 47.5 x 10^ 
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Figure 1: Radial ordering around hydrophobic particles, a Radial correlation 
function for water centers; (inset) corresponding relative density fluctuations, (bottom) 
Distribution of the number of total (dashed lines) and hydrogen-bonded (full lines) contacts 
per water molecule located in the b first and c second hydration shell with the corresponding 
bulk water values for comparison. The vertical lines denote expected values. The solute 
sizes are: r(sl)=1.4 A, r(s2)=2.1 (= 1.5 x r^^o) A, r(s3)=2.52 (= 1.8 x r^';^) Aand r(s4)=2.8 
(= 2 X r^*o) A(see Methods for details). 

Figure 2: Mutual orientational ordering of water molecules in close contact. 

Difference in the distribution of the dipolar order parameter with respect to the bulk, 
Ap{D) = p{D) — pb{D), in the first a and second b hydration shell. The dipolar order 
distribution in bulk water is shown in the inset of a. c Joint probability of finding two 
molecules at a distance r apart having a (pair) mutual order D2 (see Eq. [T] for definition, 
with N=l). 

Figure 3: Fluctuation and correlation entropy differences per tagged pair/cluster 
in the first hydration shell with respect to the bulk are given in the Tables. Cor- 
relations between differently interacting pairs of water molecules are depicted 
by colour in schematics and Tables. Only the tagged molecule necessarily lies in the 
first hydration shell. The crosses denote the tagged molecule which is used for the local- 
ization of a given cluster and the dashed black lines denote HBs. a Absolute and relative 
entropy difference of ffuctuations of HB (blue shading), non-HB (green shading) and total 
interaction energy with nearest neighbors (black frame). The slight increase of the entropy 
of total interaction energy ffuctuations accompanied with the simultaneous decrease of en- 
tropy of HB and non-HB energy ffuctuations is indicative of higher order correlations, b 
3-body correlation entropy of geminal (pairs share the tagged molecule) HB pairs (blue 
shading) and a geminal HB and non-HB pair (black frame), c 4-body correlation entropy of 
HB-bridged HB pairs (blue shading), non-HB-bridged HB pairs (black frame) and a vicinal 
(pairs do not share the tagged molecule) tagged HB and non-HB pair (magenta frame), 
d 4-body correlation entropy of vicinal HB-bridged non-HB pairs (green shading), vicinal 
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non-HB-bridged non-HB pairs (magenta frame) and a vicinal tagged non-Hb and HB pair 
(black frame). The upto 2200% increase of many-body correlations with respect to bulk 
water confirms their crucial role in mediating the hydrophobic effect. 

Figure 4: Emergence of pronounced many-body correlations and their physi- 
cal origin, a In bulk water the fluctuations are correlated predominantly up to the pair 
level. In the vicinity of hydrophobic solutes significant many-body correlations emerge in 
order to compensate for reduced fluctuations in the local electrostatic field as a result of the 

exceedingly small field of the hydrophobic solute. Different colours denote various orders 
of correlated fluctuations, b Distribution of magnitudes of the local electrostatic field ex- 
perienced by a water molecule at a given location. Black lines denote distributions in bulk 
water and symbols represent results for water molecules around apolar solute. If one nearest 
neighbor water molecule is omitted from the calculation of the local field in the bulk, this 
leads to uncompensated fluctuations shown by the red distribution. Results are given as 
relative deviations from the average bulk value \Eh\. Note that the strength of the local fleld 
experienced by water molecules in the first and second hydration shell is almost identical 
to the one in bulk water, c The coordinate frame used in the calculation of angular dis- 
tributions of local electrostatic fields. The radial projection of the water molecule position 
onto a spherical surface {9', </?') defines the origin and orientation of a secondary coordinate 
system (blue), such that the secondary z(^', </?'), x{9',ip') and y{9',(p') axes correspond to 
the local radial, polar and azimuthal directions. Angular coordinates in the secondary frame 
{9{9' ,ip'},ip{9' are used to describe the distribution of the local field over the surface 
of the unit sphere, d Distribution of orientations of the local electrostatic field in the first 
hydration shell of apolar solutes and in bulk water. In the first hydration shell the orienta- 
tional fluctuations of the local fleld are conflned to the plane containing the normal to the 
solute surface. The distribution in the bulk as well as in the second hydration shell (not 
shown) is entirely isotropic. 
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FIG. 2. 
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si s2 s3 s4 


A.S[10-2tnl 

1 Asys^i%i 


-1.2 0.4 2.6 1.9 
0.2 1.1 0.7 






ASIlO-^itfl] 
AS/SB[%] 


-2.8 -2.6 -2.2 -2.1 
-1.9 -1.7 -1.5 -1.4 


AS/S''l%] 


-4.8 -4.8 -9.1 -8.3 
-2.6 -2.6 -4.9 -4.4 







si s2 s3 34 \ 


ASJ10-"*b] 
A.S,/Sf[%] 


3.9 3.8 4,8 3.8 1 
910 884 1120 884 | 


a.s,[io-3*b1 

ASJSf[%] 


4.7 4.4 6,0 4.1 
72 68 92 63 






AS,,[10-2%| 
AS,/Sf[%] 


0.8 0.9 1.2 0.9 
145 164 218 164 






si s2 s3 s4 


|AS,[10-2tel 
1 AS,/Sf|%] 


2.2 4.3 5.6 4.1 
149 291 378 277 


ASJlO-^tfi] 

Us,/sf [%| 


0.9 1.3 1.6 1.2 
1216 1757 2162 1622 


|AS<.[10-=fc«] 
1 A.S,:/.S,"[%| 


0.7 0.9 1,2 0.9 
876 1125 1500 1126 



FIG. 3. 
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FIG. 4. Fig 4 
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I. PROBABILITY DENSITIES AND THE MAPPING ONTO POTENTIAL EN- 
ERGY SPACE 



Throughout the paper as well as this section we strictly discriminate between a random 
variable, which is denoted by uppercase letters, and a value it can take on (lowercase letters). 
The random variable transformation theorem is used to map the joint configurational prob- 
ability density, onto its functionally dependent joint probability density for pair potential 
energies. Thereby molecular positions and orientations, — {R^,<x^), are mapped onto 
a vector of interaction energies, U: 

p{u) = I du;^P(a;^)5(w - /(u;^)), (1) 

The functional relation S{u — /(w^)) contains class indicators (close contact, hydrogen- 
bonded, etc.) and the the appropriate averaging operation accounting for the indistin- 
guishability of pairs within a certain class. 

If we denote the solute position as Rg the general functional expressions 5{u — f{u)^)) 
can be written as follows. 

In the scalar case {i.e. density p{Ux)) we have: 

S{u- f{u;'')) = —^ ^ (2) 

/dw^P(aJ^) YEij^iXd{<^i,<^j)XHB{uJi,uJj)^Rr^i^,RmaA\ri - VsDj 

where Xdi^i, ^j) is the indicator function for the distance cutoff for nearest neighbors 

xM,ojj)^{l-H{\r,-r,\)), (3) 

H{x) being the Heaviside function and XHB{<^i,<^j) is the indicator function for hydrogen 
bonding (which is left out in the case of non-hydrogen bonded neighbors) 

XHsicOi, <^j) = H{d[OHi ■ ■ ■ Oj] - dmin) + H{d[OH^ ■■■Oi]- ^min), (4) 

where i^fOifj ■ ■ ■ Oj] is the angle between the OH bond of the donor molecule and vector 
connecting the proton of the donor molecule and the acceptor oxygen, 'dmin is the lower 
cutoff for the hydrogen bonding angle. Ila,b{x) represents the boxcar function 



if I I F^min ) 

1 if Rfnin — \ fi fs\ — Rmaxi 
if Vi Ve ^ Rm,ax • 



which locahzes the tagged molecule z in a certain hydration shell (first, second shell or bulk). 

In the non-scalar case p{U) is a joint distribution. Then we have, for example for a 
joint distribution for a 3-body or 4-body density (or equivalently for a pair density of pairs), 
p{Uox, Uoy) or p{Uox, Ujz) respectively: 

5(w-/(u;^)) ^ ^'^'^ ^^"^ ~ ~ ^'')xIx'Jx%BX'IiB^R^..,Rr..M^)a - SuS,iS'J,,)Q^f^ 

J du;^P(u;^) xM'x&BX'l^B^R...r.,Rr.M{l - SuW,,)Q'r} ' 

(6) 

where the explicit coordinate dependences (fij, flj) were replaced by superscripts di — 
\ri — rs\ and Sij is the Kronecker delta. The primed Kronecker delta denotes that i may or 
may not be equal to k, depending on whether we take geminal or vicinal pairs of molecules, 
that is, 3- or 4- particle correlations (for the explicit definition of geminal and vicinal pairs see 
main text). For geminal pairs the primed Kronecker delta is omitted. Q^J^'' is the connectivity 
indicator function and defines the topological relation between both pairs. Two pairs may 
be connected with a hydrogen bond or not. 

Corresponding equation for higher dimensional probability densities can be constructed 
accordingly. 



II. CALCULATION OF THE ORINENTATIONAL DISTRIBUTION OF THE LO- 
CAL ELECTROSTATIC FIELD 

The systems under investigation all have spherical symmetry around the center of the 
solute. Thus, all points in the space-fixed frame (denoted by prime) with with the same 
radius r' are positionally equivalent. As we classify the water molecules according to their 
position with respect to the solute as belonging to the first and second hydration shell or the 
bulk, we can drop the radial component when dealing with molecules in a particular shell 
and simply retain the angular variables {6',ip'), where 9' = arccos(^) and Lp' = arctan(^). 
We choose to describe the orientational distribution of the instantaneous electrostatic filed 
orientation in a given point of the space-fixed frame as a distribution over the surface of a 
unit sphere. In order to take into account the equivalence of positions at given radius in the 
space-fixed frame we set the secondary local frame in the following manner. The secondary 
X, y and z axes are chosen to coincide with the polar, azimuthal and radial directions at 
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{9',ip'). Specifically, this amounts to the following basis: 



^ cos (p' sin 6' ^ 



sin (/?' cos 6' 
-sine' 



^ — sin (f' ^ 



\ 



cos 




^ cos sin ^' ^ 



I 



sin sin Q' 
cos ^' 



The local electrostatic field at a point r' = (x', y', z') due to a collection of point charges 
Qi with coordinates r" is -E(r*') = ^f=i ^rpz^pm-- The unit Cartesian components of the field 



in the secondary frame at (r', 9', (p') are thus: 



[E{r') ■ y) 
\E{r')\ 



[E{r') ■ z) 
\E(r')\ 



After a trivial transformation to spherical coordinates in the secondary frame, we may write 
for the orientational distribution of the electrostatic field experienced by a water molecule 
in shell di (first, second or bulk): 



p{9,(p)[dj 



Rr, 



■mm j-'t-max 



III. GIBBS-SHANON AND CORRELATION ENTROPY DIFFERENCES 



Wc denote a tagged molecule with 0, its neighbors with i, j, k; i ^ j ^ k ^ . . ., and the 
nearest neighbors of the nearest neighbors with primes. A 3-body correlation entropy of a 
tagged molecule and its two nearest HB neighbors is thus, for example, S^^^J^q, and the 
4-body correlation entropy of two adjacent HB pairs is S^^-'^jj^, if j is not hydrogen-bonded 
to 0, and S'^^b-Vb' J ^ HB neighbor of 0. . 
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TABLE I. Entropy differences per water molecule in the first and second shell for various solute 
sizes. ASi stands for 5* — 5^ and the asterisk denotes that molecule j is hydrogen-bonded to the 
tagged molecule 0. 





sl 1** 


sl 2"^^ 


s2 1** 


s2 2"'^ 


s3 1** 


s3 2'"^ 


s4 P* 


s4 2^^°' 




-7.9 


-6.4 


4.5 


2.0 


-4.6 


-6.1 


5.6 


3.3 


asd/sE[io-'%] 


8.6 


6.9 


-4.9 


-2.2 


5.0 


6.7 


-6.2 


-3.6 


AStot^c[W~^kB] 


-1.2 


-0.3 


0.4 


0.2 


2.6 


-0.3 


1.9 


0.2 


^Stot-c/ Stot-c\%] 


-0.5 


-0.1 


0.2 


0.1 


1.1 


-0.1 


0.7 


0.1 


A5i/B[10-2/cB] 


-2.8 


-1.1 


-2.6 


0.4 


-2.2 


0.07 


-2.1 


-0.5 


AShb/S^M 


-1.9 


-0.7 


-1.7 


0.3 


-1.5 


0.05 


-1.4 


0.3 




-4.8 


-0.5 


-4.8 


-0.1 


-9.1 


-0.3 


-8.3 


-0.9 


ASnHB/ S^glVo] 


-2.6 


-0.3 


-2.6 


-0.05 


-4.9 


-0.16 


-4.4 


-0.5 


^'^HB-HB/ '^HB-HBi^^l 


5.1 
94 


1.7 
31 


4.8 
91 


1.3 

25 


6.7 

126 


1.9 

36 


5.2 

98 


1.8 

34 


^S'^B'-^nHBi^^-'kB] 

AqiO-Oj ,qiO-Oj;B rqy-i 
^'^HB-nHB/ '^HB-nHB\-^°\ 


6.9 

812 


2.0 

235 


6.0 

706 


1.8 
212 


7.4 

871 


2.2 

256 


6.1 

718 


1.7 

200 


AS'^^'U^O-^ks] 

A qiO-jk' 1 qiO-jk;B r(v i 
^'^HB-HBI •^HB-HBV^^\ 


3.9 
910 


1.2 

440 


3.8 

884 


1.1 

260 


4.8 
1120 


1.3 

320 


3.8 

884 


1.1 

256 


^S'^7-HB[^0-^kB] 

AqiO-jk'* / qiO-jk*;Bro^-i 
^^HB-HbI ^HB-HB [/OJ 


4.7 
72 


1.2 

18 


4.4 

68 


0.9 

14 


6.0 

92 


1.4 
22 


4.1 

63 


1.4 

22 


^S^ntuHBi^^-'kB] 

A qiO-jk' 1 qiO-jk;B rcr/-! 
'-^'^nHB-nHBI '^nHB-nHBy^°\ 


2.2 

149 


0.6 

41 


4.3 
291 


1.1 

74 


5.6 

378 


1.3 

88 


4.1 

277 


0.9 
61 


^"^nHB-nHBy^^ '^B\ 

xqiO-jk'* /qiO-jk*;B r(wi 

'-^'^nHB-nHBl '^nHB-nHBV^\ 


0.9 
1216 


0.3 
405 


1.3 

1757 


0.3 

405 


1.6 

2162 


0.4 

541 


1.2 

1622 


0.3 

405 


^^H-J-nHB[^^-'kB] 

\qiO-jk' ,qiO-jk;B r(yi 
^'^HB-nHBl '-'HB-nHBl^'^i 


0.7 

875 


0.2 

250 


0.9 

1125 


0.2 

250 


1.2 

1500 


0.3 

375 


0.9 

1125 


0.2 

250 


^qiO-jk'* ,qiO-jk*;B rfw-i 

'-^'^HB-nHB 1 '^HB-nHB L 


0.8 
145 


0.2 
26 


0.9 
164 


0.2 
26 


1.2 
218 


0.3 
55 


0.9 
164 


0.3 
55 
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